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A Polynomial Time Algorithm to Compute an Approximate 

Weighted Shortest Path 

R. Inkuhr Sanjiv Kapoor* * 


Abstract 

We devise a polynomial-time approximation scheme for the classical geometric problem of 
finding an approximate short path amid weighted regions. In this problem, a triangulated region 
V comprising of n vertices, a positive weight associated with each triangle, and two points s and 
t that belong to V are given as the input. The objective is to find a path whose cost is at most 
(1 + e)OPT where OPT is the cost of an optimal path between s and t. Our algorithm initiates 
a discretized-Dijkstra wavefront from source s and progresses the wavefront till it strikes t. This 
result is about a cubic factor (in n) improvement over the Mitchell and Papadimitriou ’91 result 
[8], which is the only known polynomial time algorithm for this problem to date. Further, with 
polynomial time preprocessing of V, a map is computed which allows answering approximate 
weighted shortest path queries in polynomial time. 


1 Introduction 

The problem of computing a shortest path in polygonal subdivisions is important and well-studied 
due to its applications in geographic information systems, VLSI design, robot motion planning, etc. 
A survey of various shortest path problems and algorithms may be found in Mitchell [9]. In this 
paper, we devise an algorithm for the weighted shortest path problem |8j: given a triangulation V 
with 0(n) faces, each face associated with a positive weight, find a path between two input points 
s and t (both belonging to V) so that the path has minimum cost among all possible paths joining 
s and t that he on V. The cost of any path P is the sum of costs of all line segments in P , whereas 
the cost of a line segment is its Euclidean length multiplied by the weight of the face on which it 
lies. The weighted shortest path problem helps in modeling region specific constraints in planning 
motion. 

To compare the time complexities of various algorithms in the literature, we use the following 
notation - n: number of vertices defining V; l m ax■ length of the longest edge bounding any face 
of V] N : maximum coordinate value used in describing V ; w m ax■ maximum non-infinite weight 
associated with any triangle; w m i n : minimum weight associated with any triangle; 6 m i n : minimum 
among the internal face angles of V; and, p: ratio of w m ax to w m i n . (Note that the same notation 
is used in later parts of the paper as well.) 

Mitchell and Papadimitriou [8] presented an algorithm that finds an approximate weighted 
shortest path in 0(n 8 lg hKh) time. Their algorithm essentially builds a shortest path map for s by 
progressing continuous-Dijkstra wavefront in V using the Snell’s laws of refraction. By introducing 
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m equi-spaced Steiner points on each edge of V and building a graph spanner over these points, 
Mata and Mitchell [7] devised a preprocessing algorithm to construct a graph spanner in 0{kn 3 ) 
time, where k = O(jg-), on which (1 + e)-approximate weighted shortest path queries are per¬ 
formed. Lanthier et al. |B] independently devised an 0(n' 5 ) time approximation algorithm with an 
additive error of 0(l max w max ) by choosing m = n 2 . Instead of uniform discretization (as in 0 ), 
Aleksandrov et al. [2113] used logarithmic discretization and devised an lg(^y-) lg j lg 

time approximation algorithm. Sun and Reif m provided an approximation algorithm, popularly 
known as BUSHWHACK, with time complexity 0(p^ lg (Nfj,) lg f lg ^). Their algorithm dynam¬ 
ically maintains for each Steiner point v, a small set of incident edges of v that may contribute 
to an approximate weighted shortest path from s to t. More recently, Cheng et al. [3] devised 
an approximation algorithm that takes ( 'J[ kn+k teR'A) }g 2 PH ) time, where k is the smallest integer 
such that the sum of the k smallest angles in V is at least ir. The query version of this problem is 
addressed in ! 2 i is uni Euan]. Further, algorithms in Cheng et al. [5] handle the case of measuring 
the cost of path length in each face with an asymmetric convex distance function. 

Our contribution 

The time complexities of each of the above mentioned solutions, except for [8], are polynomial in n 
as well as in parameters such as e, 6min, N and ji. Hence, strictly speaking, these algorithms are not 
polynomial. Like [8], this paper devises an algorithm that is polynomial in time complexity. The 
time complexity of our algorithm is 0(n 5 ( lg (r)(lg ^)), which is about a cubic factor improvement 
from [8j. As established in pB], there are kl(n 4 ) events that need to be handled in order to find the 
interaction of shortest path map with the V. Our algorithm takes first steps to provide a solution 
that is sub-quadratic in the number of events. 

This result uses several of the characterizations from |8j in order to design a simple and more 
efficient algorithm. Every ray is a simple path in V. Our approach discretizes the wavefront in 
the continuous-Dijkstra’s approach by a set S of rays whose origin is source s. These rays are 
distributed uniformly around s. As the discrete wavefront propagates, a subset of the rays in S 
are progressed (traced) further while following the Snell’s laws of refraction. Each of these subsets 
of rays is guided by two extreme rays from that subset, i.e., all the rays in that subset lie between 
these two special rays. Essentially, each such subset represents a section of the wavefront and is 
called a bundle. For any vertex v in V, whenever such a section of the wavefront strikes a vertex 
w, we initiate another discrete-wavefront (set of rays) from w. We continue doing this until the 
wavefront (approximately) strikes t. In summary our contributions are: 

1. A discretized approach to propagating wavefronts using bundles. 

2. An algorithm that computes an (1 + e) -approximate weighted shortest path from s to t in 
0(n 5 (lgf)(lg^)) time. 

3. Further, we preprocess V in 0(n 5 (lg ”)(lg ^=)(lg IV)) time to compute a data structure for an¬ 
swering single-source approximate weighted shortest path queries in 0(n 4 (lg ")(lg ^=)(lgIV)) 
time. 

Section [2] lists relevant propositions from [8] and [TO] . and defines terminology required to 
describe the algorithm. Section [3] outlines the algorithm while introducing few structures used in 
the algorithm. In Section]!] we bound the number of rays. The details of the algorithm are provided 
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in Section [5j Section [6] describes an interpolation scheme to improve the time complexity of the 
algorithm. Section [7] argues for the correctness and analyzes the time complexity of the algorithm. 
And, the conclusions are given in Section [8j 

2 Preliminaries 

We define the problem using the terminology from [8j. 

We assume a planar subdivision V, that is polygonal and specified by triangular faces. Each 
face / has a weight Wf associated with it. We denote the weight of a face / (resp. edge e) with Wf 
(resp. w e ). For an edge e shared by faces f and f", the weight of e is defined as min (wf,Wfn). 

A path is a continuous image of an interval, say [0,1], in the plane. A geodesic path is a path 
that is locally optimal and cannot, therefore, be shortened by slight perturbations. An optimal 
path is a geodesic path that is globally optimal. The general form of a weighted geodesic path is a 
simple (that is, not self-intersecting) piecewise linear path that goes through zero or more vertices 
while possibly crossing a zero or more edges. The Euclidean length of a line segment l is denoted 
by ||Z||. Let P be a geodesic path with line segments ■■■ such that l t lies on face /*, for 
every i in [1, k]\ then the weighted Euclidean distance (also termed the cost ) of path P is defined as 
the EJI*i||un¬ 
weighted Shortest Path problem: Given a finite triangulation V in the plane with two points 
s (source) and t (destination) located on V, an assignment of positive integral weights to faces of 
V, and an error tolerance e, the approximate weighted shortest path problem is to determine a path 
P from s to t that lies on V such that the cost of P is at most (1 + e) times the cost of an optimal 
path from s to t. 

We assume that all input parameters of the problem are specified by integers. In particular, all 
vertices have non-negative integer coordinates. Further, we assume that s and t are two vertices 
oiV. A ray in our algorithm is a piecewise linear simple path in V such that endpoints of each of 
the line segments that it contains lie on the edges of V. Let f be the first face traversed by a ray 
r' and let f" be the first face traversed by a ray r". The angle between r' and r" is defined as the 
angle between the vectors induced by f n r' and f" n r". Unless specified otherwise, all angles are 
acute. Let it, v\ and it) be three unit vectors that originate from a point a in M 2 such that the 
vectors v\_, v) he in the same half-plane defined by the line induced by 1 1. The cone C(a, ut, v)) is 
the set comprising of all the points that are positive linear combinations of Fi and v). 

A few facts from the literature 

First we list few propositions, definitions and descriptions from Mitchell and Papadimitriou [Sj 
tailored for our purpose. 

A sequence of edge-adjacent faces is a list, (/i,/ 2 , • • • ,/fc+i), of two or more faces such that, 
for every i, face fi shares edge e* with face fi+i- Further, the corresponding sequence of edges 
£ = (ei,e 2 , • • • ,ek) is referred as an edge sequence. When a geodesic path P crosses edges in £ in 
the order specified by £ and without passing through any vertex, then £ is the edge sequence of 
path P. A geodesic path po,pi,... ,Pk in V, is termed a ray r as it behaves similar to a ray of light. 
Further, when the geodesic path is extended to add a point pu+i in V , by the line segment PkPk+i , 
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the ray r is said to be traced (or, progressed) to Pk+i- Any point or line segment in V from which 
at least one ray is initiated is termed a source. 

Let f and f" be two faces with shared edge 
e(v',v"). (See Fig. [lj) Let l',1" be two suc¬ 
cessive line segments along a ray r with l' lying 
on face f and l" lying on face f" with point 
y 6 l' fl l" n e. Let v be a vector entering y and 
containing l' and let i/' be a vector with origin 
y and containing l". Also, let it be a vector 
normal to edge e, passing through point y to some point in face f". The angle O’ between it and 
v’ is known as the angle of incidence of r onto e. And, the angle 6" between it and i/' is known 
as the angle of refraction of r from e. 

When wp > Wf", the critical angle of e, 
denoted by O c {f,f") is sin _1 (^). (When f 

and f" are understood from the context, the 
critical angle of e is also denoted by 6 C .) (See 
Fig. [2]) If O’ < O c (f',f"), then the angles O’ 
and 6" are related by Snell’s law of refraction 
with wpsinO' = Wf" sva.6". Since there does 
not exist a geodesic weighted shortest path with 
O' > O c (f',f"), we only need to consider the case in which O' equals to 6 c (f,f"). Let y' be the 
point of incidence of ray r' from face f onto e with angle of incidence 0 c (f , f"). Also, let it be a 
vector normal to edge e, passing through point y' to some point in face f". Then the geodesic path 
travels along y'v" for some positive distance before exiting edge e back into face f, say at a point 
y" located in the interior of y'v" , while making an angle — 9 c (f , f") with it. We say that the path 
P is critically reflected by edge e and the line segment y'y" is termed a critical segment of path P 
on e. Sometimes, we also say y'v" is a critical segment corresponding to the critcial incidence of a 
ray at y' from face f. The point y' (the closer of the two points {y', y"} to s) is known as a critical 
point of entry of path P to edge e and y" is known as the corresponding critical point of exit of 
path P from edge e. 

The following propositions from the literature are useful for our algorithm. 

Proposition 1 (Lemma 3.7, [8]I Let P be a geodesic path. Then either (i) between any two 
consecutive vertices on P, there is at most one critical point of entry to an edge e, and at most one 
critical point of exit from an edge e' (possibly equal to e); or (ii) the path P can be modified in such 
a way that case (i) holds without altering the length of the path. 

Let f and f" be two faces with shared edge e. For any point x G e, a locally f -free path strikes 
x from the exterior of face f and is locally optimal. 

Proposition 2 (Lemma 7.1, [S]) For a face f of V , let P be a shortest locally f-free path. Let 
P' be a sub-path of P such that P' goes through no vertices or critical points. Then, P' can cross 
an edge e at most 0(n) times. Thus, in particular, the cardinality of any edge sequence of path P 
is 0(n 2 ). 

Proposition 3 (Lemma 7.4, [8]) There are at most 0(n ) critical points of entry on any given 
edge e. 
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Figure 2: Illustrating critical reflection 
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Finally we have a proposition that provides us the non-crossing property of (weighted) shortest 
paths. 

Proposition 4 (Lemma 1, jlOj ) Any two geodesic shortest paths that originate from the same 
point in V cannot intersect in the interior of any weighted region of V. 

3 Algorithm outline 

In our algorithm, we progress a discretized-Dijkstra wavefront, as the traditional approach of pro¬ 
gressing a continuous-Dijkstra wavefront results in a complicated algorithm, This wavefront is de¬ 
fined and expanded using rays. Every ray is a geodesic path in V from its point of origin to a point 
on the wavefront. The wavefront that defines the locus of points at weighted Euclidean distance d 
from the source is identified by points on the the rays that are at weighted Euclidean distance d 
from the source. Note that geodesic paths in V are piece-wise linear. We initiate a set TZ(s) of rays 
whose origin is s. The rays in 7 Z(s) are ordered according to their counterclockwise angle with the 
positive x-axis and are uniformly distributed around s. These rays together describe a (discrete) 
wavefront initiated at s. Let r' and r" be two rays in IZ(s). Two rays are termed successive when 
they are adjacent in the ordering of rays around any point from which those rays are initiated. The 
number of rays in TZ(s) is defined by the angle 6 between successive rays; by expressing the value 
of 5 in terms of e later (in Section [4]), our algorithm ensures an (1 + e)-approximation. When the 
wavefront strikes any vertex v, an ordered set lZ(v) of rays are initiated from v, unless vertex v is 
the destination t itself. 

3.1 Types of rays 

We next discuss the types of rays that form part of our wavefront. Only a subset of the rays 
that have been initiated are used in the propagation of the wavefront; a ray that is considered for 
propagation is said to have been traced. We propagate the rays as described below: let e = (v 1 , v") 
be a common edge between faces f and f". Consider a ray r G 7£(u) that is traced along face f 
and suppose it strikes e at a point y' G e. 

If the angle of incidence of r onto e is less than 0 c {f, f"), then r refracts onto face f" with the 
angle of refraction defined by Snell’s law of refraction. The case in which the angle of incidence of 
r onto e is greater than 8 c (f, f") occurs only when wy > wyr. In this case, we do not propagate 
r further as r would not be part of a weighted shortest path to t. 

If the angle of incidence of r onto e is equal 
t° Oc(f, f")i then a weighted shortest path that 
uses the ray r propagates along e for a pos¬ 
itive distance before critically reflecting back 
into face f itself. Let k = y'v" be the critical 
segment corresponding to this critical incidence 
of r onto e where y' is the point of incidence of 
r on e. Since a weighted shortest path can be 
reflected back from any point on the critical segment k of e, our algorithm initiates rays from a 
discrete set of evenly spaced points on k. The number and position of points from which these rays 
are generated is again a function of e. Let it be a vector normal to edge e, passing through point 
y' to some point in face f". These rays critically reflect from e back into face f while making an 
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angle —9 c (f',f") with if. (See Fig. [3|) We let 1Z(k) be the ordered set of rays that originate from 
k. and are ordered by the distance from y'. Two rays r' and r" in 7Z{n) are successive whenever 
there is no ray in 1Z(k) that occurs between r' and r" in the linear ordering along k. 

We next introduce an additional category of 
rays. To account for the divergence of rays and 
to ensure that the wavefront is adequately rep¬ 
resented, we further initiate an angle ordered 
set of Steiner rays , denoted by lZ(y '), from each 
critical point of entry y' onto face f". (See 
Fig. |4j) Hence, a critical point of entry is 
also termed as a critical source. These rays are 
motivated as follows: a pair of successive rays 
when traced along an edge sequence can diverge 
non-uniformly, with the divergence being large 
especially at angles close to critical angles. To 
establish an approximation bound, instead of having a set of rays that is more dense from every 
source, we fill the gaps by generating Steiner rays from critical points of entries. This helps in 
reducing the overall time complexity. Let v' and v" be the endpoints of e. For any source u, let 
r, r' be successive rays in 7 Z(u) such that r is refracted and r' is critically reflected. Let y (resp. 
y') be the point at which r (resp. r') is incident to e = ( v',v Let O' be the angle of refraction 
of r. Also, let V 2 be a vector with origin at y' and that makes an angle 9' with if. As a result of 
discretization, there may not exist rays in 1Z(u) that intersect the cone C(y', y'v\ Hence, to 
account for this region devoid of rays, an ordered set 7 Z(y') of rays are initiated from y' that lie 
in the cone C(y\ y'v\ V 2 ), ordered by the angle each ray makes with respect to if. Every ray in 
7 Z(y') is termed a Steiner ray. 

To recapitulate, there are three sets of rays that are used to define the discrete wavefront: 

• {7 Z(v) | v eV} 

• {7Z(n) I k is a critical segment obtained when a ray r £ 7Z(v) strikes an edge e at a critical 
angle} 

• {7 Z(y') | y' is a critical point of entry for a ray r £ 7 Z(v) for some v £ V} 

3.2 Ray bundles 

We approximate the expansion of the continuous-Dijkstra wavefront by tracing rays in pairs (angle 
between the rays being less than ir) so that each pair represents a section of the continuous- 
Dijkstra wavefront that lies between the pair at corresponding distance d along the chosen pair of 
rays. Appropriately chosen pairs suffice to represent the behavior of this section. The rays could 
be initiated from either of these: (i) a vertex v, (ii) a critical source c, or (iii) a critical segment k. 
Let o be one such source of rays. We partition the set of rays with origin o as follows: Let B be a 
maximal set of successive rays in 7 Z(o) such that all rays in B cross the same edge sequence £ when 
traced from the source to the current state of the discrete wavefront. Then B is said to be a ray 
bundle of 7Z(o). Furthermore, £ is the edge sequence associated with B. Let r / and r" be two rays 
in B intersecting an edge e £ £ such that the line segment defined by the points of intersection 
of r' and r" with e intersects every other ray from B (when traced). Then the rays r' and r" are 



Figure 4: Illustrating rays initiated from a critical source 
y' and the corresponding sibling pair r, r 3 
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extremal rays of the bundle and are termed as the sibling pair of ray bundle B. For any ray bundle 
B, instead of tracing all the rays in B, we trace only the sibling pair of B. This helps in reducing 
the number of event points, hence the time complexity. 

The ray bundles and corresponding sibling 
pairs are updated as the rays progress further. 

For the sibling pair r', r" of ray bundle B , when 
the sequence of edges intersected by r' changes 
from the sequence of edges intersected by r" for 
the first time, we need to split the bundle of 
rays into two and determine new sibling pairs. Figure 5: Split of a ray bundle 

A binary search among the rays in B is used to trace the appropriate rays across the edge sequence 
of B and form new sibling pairs. Let abc be a face / in V. (See Fig. 0) Suppose rays r' and r" 
are siblings in a bundle B before they strike edge ab of /. However, when they are traced further, 
suppose r' strikes edge ac and r" strikes edge be. At this instance of wavefront progression, the 
sequence of edges that are intersected by r' and r" differ for the first time, and hence r' and r" do 
not belong to the same bundle from there on. Using binary search over the rays in B, a pair of 
successive rays, say and r”, are found to determine new sibling pairs: r',r[ and r ", r'{. Further, 

the ray bundle B is split accordingly. And, the wavefront is said to strike c. The details of the 
algorithm to find rays r\ and r'[ is described later. 

The sets of rays in TZ{y') and 77 (k), where y' is a critical point of entry and n is a critical 
segment, are handled similarly. The rays in 77 (k) have the characteristic that the paths of any two 
rays in IZ(n) are parallel between any two successive edges in the edge sequence associated with 
the corresponding ray bundle. 

3.3 Tree of rays 

Let v be a vertex in V. And let TZ(v) be the set of rays initiated from v. A ray may lead to the 
initiation of critical sources as described before. Let Si be the set of critical sources such that for 
any source q in Si, Steiner rays are initiated from q when the discrete wavefront that originates at 
v strikes q. Further, let S 2 ,..., S,_ 1 , Si, ..., S& be the sets of critical sources such that for every 
q' E Si and 2 < i < k, critical source q' is initiated due to the critical incidence of a Steiner ray 
from a source in Si-±. We organize the sources in the set S(v) = {u} U (U jSj) into a tree of rays , 
denoted by Tr(v). Each node of Tr(v) corresponds to a source in S(v): more specifically, v is the 
root node and every other node in Tr{v) is a distinct critical source from U jSj. For any two nodes 
u E Si-i, w E Si in Tr(v), u is the parent of w in Tr(v) if and only if the critical source w, located 
on an edge e, is initiated when a ray from u strikes e at the critical angle for edge e (and the 
point of critical incidence is w). The in-order traversal of the tree Tr{v) provides a natural order 
on the set of rays. A subset S of rays are termed successive whenever rays in S are a contiguous 
subsequence of the ordered set of rays 

3.4 Ray Bundles and Sibling Pairs 

In order to partition the rays in lZ(v), where v is a vertex, or in 77(/v), where n is a critical segment, 
we generalize the definitions of sibling pair and ray bundles. 

We consider the rays in IZ(v). Let £ = {e±, e 2 , • • •, e^} be an edge sequence. For any j E [1, k\, 
we say ej, ..., e* is a suffix of edge sequence £. Two rays belonging to Tr{v) are siblings whenever 
the edge sequence associated with one of them is a suffix of the edge sequence of the other. A 
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maximal set S of successive rays initiated from the nodes of 7r(v ) that are siblings to each other is 
a ray bundle. For two rays r', r" belonging to a ray bundle B and points p' = r’ n &k,p" = r" n e^, 
the rays r' and r" are termed as the sibling pair of B whenever the line segment p'p" intersects 
every ray in B. A similar definition holds for rays in 1Z(k). 

Therefore, there are two kinds of sibling pairs possible: Either (i) both the rays in a sibling pair 
originate from the same source or the sources of both the rays in a sibling pair belong to a tree of 
rays, Tr(v ) or (iii) both the rays in a sibling pair originate from the same critical segment. A ray 
bundle of the first category is denoted as a ray bundle of 7r(v). 

Let r',r" be a sibling pair of rays in a bundle B in Tr(v). Let u,u> G Tr(v) be the respective 
origins of r' and r". Also, let v' be the least common ancestor of u and w in Tr{v). We define v' to 
be the root of B. The path P in tree 7r(v) from v' to u comprising of a sequence of critical sources 
is the critical ancestor path of r' with respect to r". Note that the critical ancestor path of a ray is 
always defined with respect to a sibling pair. Given that rays in the bundle B are successive rays, 
nodes on the critical ancestor paths from u and w to v 1 in the subtree Tr(v) provide all the rays in 
the bundle B. 

We next show a property of the rays in a bundle. We define two rays to be pairwise divergent 
during traversal of an edge sequence £ if they do not intersect each other during the traversal. A 
bundle of rays is divergent during traversal of £ if all rays in the bundle are pairwise divergent 
during the traversal. 

Lemma 3.1 Let B be a ray bundle of Tr{v) and let £ be the edge sequence associated with rays in 
B. The rays in B are pairwise divergent whey they traverse across the edge sequence £. 

Proof: Let v' be the root of a bundle. The rays in B that start from v 1 , termed R(v'), clearly 
diverge during the traversal of the sequence of edges £, since divergence is preserved on refraction. 
Additionally, consider rays in B that are generated from a critical source, w, termed R(w). Each 
pair of rays, a steiner ray r originating at w and a ray in R(v'), are pairwise divergent, by con¬ 
struction. Similarly, every pair of ray in R(w) is pairwise divergent. These rays remain pairwise 
divergent after refraction during propagation across edges in the sequence £. In general, let R be 
a bundle of rays that is divergent during traversal of a subsequence of £. These rays are generated 
from a vertex and/or a set of critical sources. Suppose the rays in the bundle strike an edge ej E £ 
and refract, and also generate a critical source w' on ej. Then R U R(w'), where R(w') is the set 
of Steiner rays generated from uf, is a set of rays with every pair of rays being pairwise divergent 
during traversal of £ since every steiner ray in R(w') is pairwise divergent from each ray in R. Thus 
these rays will not intersect each other as they strike edges ej + \.. . e*,. □ 

The consequence of this Lemma is that rays in a bundle B can be angle ordered. As described 
later, this ordering and the critical ancestor path of a ray will be useful in efficiently splitting a 
sibling pair of Tr{v). 

4 Bounding the number of initiated and traced rays 

Before detailing the algorithm, we bound the number of initiated and traced rays in V. An ap¬ 
proximate weighted shortest path from source s to f can be split into sub-paths such that each 
such sub-path goes between a pair of vertices. Further, each of these sub-paths are classified into 
following types: 



• a Type-I path from one vertex to another vertex that does not use any critical segment in 
between, and 

• a Type-II path from one vertex to another vertex that critically reflects and uses at least one 
critical segment in-between. 


We next bound the approximation achieved for both of these types of paths by establishing 
properties regarding the rays. We show that given a Type-I (resp. Type-II) weighted shortest path 
P from a vertex, say u, to another vertex, say v, there exists a Type-I (resp. Type-II) path P' 
using only rays in our discretization such that P' closely approximates P. We first consider the 
divergence of successive rays that are within the same bundle. This divergence is due to refraction 
of rays. 


Lemma 4.1 Let r',r" he two successive rays in a bundle in Tr(u). Also, let q' and q" he points 
on rays r' and r" respectively, at weighted Euclidean distance d from s and lying on the same face. 
If the angle between rays r' and r" is upper bounded by ^ l -(e') n e', then the weighted Euclidean 
distance of the line segment joining q' and q" is upper bounded by de'. Here, e' is a small constant 
(expressed in terms of input parameters). 



Proof: We first show a bound on the di¬ 
vergence between a pair of successive rays in 
Tr(u ) . Let the sibling pair of rays r', r" in Tr{u ) 
traverse across the (same sequence of) faces 
/i, / 2 ... f p+ i with weights w 1 ,w 2 ,..., iv p+ i re¬ 
spectively. Let 5 be the angle between r' and 
r". Let ei, e 2 , ■ ■ ■ e p be the edge sequence as¬ 
sociated with this sibling pair, with being 
the edge adjacent to faces /* and fi + \ for ev¬ 
ery i < p. Let 0[, 9 2 ,.. ■, Op be the angles at 
which the ray r' is incident on edges ei,e 2 ,..., 
be the angles at which the ray r' refracts at edges e\,e 2 ,...,e p respectively. Similarly, let 
9[ + <5j, 0 ' 2 + 6 ' 2 , ■ ■ ■, 0' p + 5 p be the angles at which the ray r" is incident on edges e\,e 2 ,... ,e p 
respectively. And, let 0 2 + 5 2 , O 3 + £ 3 ,..., 9 P+ \ + <5 p+ i be the angles at which the ray r" refracts 
at edges ei, e 2 ,..., e p respectively. (See Fig. 0 ) For every i, if a critical angle exists for edge ei, 
assume that both 9\ and 0\ + S[ are less than that critical angle. Further, for every i, assume that 6 i 
is a small positive number less than e' (we will justify this assumption later). We will first provide 
a bound on 5 P in terms of Sj. 


Figure 6 : Illustrating the angle of refractions 
of rays r' and r" 

respectively. (See Fig. @) Let 9 - 2 ,9s ,..., 9 P -|_i 


For every integer i 6 [l,p], using Snell’s law, 

iVi sin 9\ = Wi + 1 sin 0 i+\ 


( 1 ) 


and 


Wi sin ( 9\ + b'f) = w i+ i sin (9 i+ 1 + <5j+i). 


Since = Si, 


Wi sin (#' + Si) = w i+ 1 sin (9 i+ 1 + <5j+i). 


9 



Both 5i and <5j + i are assumed to be very small; thus we approximate sin Si. cos sin dj + j. and 
cos S i+ i with the first term of the series expansion (an analysis with higher order terms reveals no 
additional benefit) and thus: 


Wi sin Q'i + wA cos 9[ = w i+ \ sin 9 i+ 1 + w i+ id i+ i cos 9 i+ i 
=> Wi5i cos 9\ = w i+ i5i + i cos 9 i+ i (/romffl 


<5j+1 — 


Wi 


COS 9': 


Letting maxi<j< p 


cos O'- 


cos Oi 


+1 


w i+ i cos 9 i+ 1 
= /?, the above leads to 

W\ 


^p+l 


VJp+1 


P p 5 1, P> 1 


(2) 

(3) 


Let 9 l c be the critical angle corresponding to faces ft and fi+i- If 9\ is close to the critical angle 
9\ then /3 grows unbounded. We thus restrict 9\ to 9 l c — Ke', where K is a constant defined as 
max i( m “siuV )• And, we note that K is upper bounded by p?. This ensures that every refracted 
ray is at an angle less than f — e 7 with respect to the normal at edge e*. Thus B < — A— n < 4. 

Z COS^ 2 ^ ) € 

The error introduced due to this assumption will be bounded later. 



Figure 7: Illustrating the construction in proving Lemma 


4.1 


(Except for d e , all the distances’ shown in the figure are 


weighted Euclidean distances.) 


To compute the weighted Euclidean distance between two points at weighted Eulcidean distance 
d from the source, first we restrict our attention to two points that are on successive rays in the 
same bundle at weighted Eulcidean distance d from the source and that lie on the same face. Let 
the two points, q' and q", both on the same face fi+u be such that they (i) lie on a pair of successive 
rays r' and r" in the same bundle in Tr(u), and (ii) are at equal weighted Euclidean distance from 
u. Let pi be the point of incidence of r' on an edge e* (which is common to f t and fi+i) with an 
angle of incidence 9[ and let it refract from p\ at an angle of refraction Further, let the ray 

r" be incident on e* at p' 2 with an angle of incidence 9[ + 5[ and refract from p' 2 at an angle of 
refraction 9 l+ \ + dj+i. (See Fig. [T}) W.l.o.g., assume that the weighted Euclidean distance from 
vertex u to p 2 is larger than the distance from u to p\. Consider an additional point p 2 located on 
ray r" such that both p\ and p -2 are in the face /,-, and the weighted Euclidean distance from u to 
either of these points is dl. Let d\ be the weighted Euclidean distance between p\ and p 2 - Also, 
points q ', q" are located on rays r 7 , r" respectively such that q ', q" are in the region fi + \ with weight 
|-i. Further, we let d be the weighted Euclidean distance from u to either of these points and let 
d 2 be the weighted Euclidean distance between q' and q". 
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We establish a bound on d 2 for i > 0. For every j, we let k ? - = ^f0^ l 5\. We first show that if 

J Wj 

d\ < 2d'/vi then d 2 < 2dKj + i. Clearly, di < 2d! K\, i.e., when i = 1. 

We first consider the case when points q' and q" are on the same face fi+\- Let p\ be the point 
that is at weighted Euclidean distance d[ from pi. Let d e be the Euclidean distance between p\ 
and p' 2 along edge e. Also, let d 2 be the weighted Euclidean distance between p[ and p 2 . 

Now by triangle inequality 

de < — (di + d[) 

Wi 

By assumption d\ < 2d'Ki and, by triangle inequality 

d2 < ^*+i(—-—I- {d\ + d'i)) 

W i+1 Wi 

< ^ + — ( 2 ^ + 4 ) 

Wi 

< 2 (d + d0)Ki -|_i 

The last inequality follows from the assumption that Kj+i > Ki > 1. 

Next consider the weighted Eulcidean distance between q' and q" , where q' is at weighted 
Euclidean distance d' 2 from p\. and q" is at weighted Euclidean distance d' 2 from p 2 . 

d'2 0: d-2 + 5i+id 2 

using the small angle approximation of distance along a circular arc subtending a small angle. 
Thus, since <5j + i < Ki + 1 , 

d 2 < 2 (d' + d[ + d' 2 )ni + 1 = 2 dn i+ i 

Since Kj = t00^ 1 5i, 0 < \ and j < n 2 , if we choose <5i < - then the weighted Euclidean 

J lU j c Z jJj 

distance d 2 is bounded by de'. □ 


The following corollary to the above eliminates the restriction for points q' and q" being on the 
same face. 


Corollary 1 Let r',r" be two successive rays in a bundle B in Tr(u). Also, let q' and q" be two 
points on r' and r" respectively, at weighted Euclidean distance d from s. If the angle between rays 
r' and r" is upper bounded by 2/j (e') n e', then the weighted Euclidean distance of the line segment 
joining q' and q" is upper bounded by de'. Here, e' is a constant (expressed in terms of input 
parameters). 


Proof: Let 5 be the angle between r' and r". Consider the locus of points at weighted Eu¬ 
clidean distance d from the root of Tr{u ) that he in the region bounded by rays r' and r". Let 
q' = Po,Pi,P 2 , ■ ■ ■ Pki o" = Pk+i be the sequence S of points of intersections of this locus with the 
edges of V. Consider the rays in the bundle B that strike points in S. Further, for 0 < i < k, let 6 i 
be the angle between rays that strike points pi and Pi+i of S. Then the Lemma 4.1 shows that the 
weighted Euclidean distance between pi and Pi+\ is 2c ^ for every 0 < i < k. The overall bound 
follows by summing up the weighted Euclidean distances’ between successive points in S. □ 


Before we proceed further, for simplicity, we establish another restriction on the angle that 
two successive rays in 7r(u ) make at their source u. To ensure that each triangle in our planar 
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subdivision gets intersected by at least two rays, we further restrict angle between any two successive 
rays so that del < Here, d = nw ma xlmax is the largest weighted Euclidean distance in V; l m ax 
and Imin are respectively the lengths of edges with maximum and minimum Euclidean lengths in 
V. Thus this establishes another condition on e', i.e. 


We can now prove the following. 


e'< 


lr 


4 nw max lmax 


(4) 


Lemma 4.2 Let u be a vertex and let v be a point on an edge e in V such that the weighted 
Euclidean distance from u to v via a Type-1 ray in bundle B is d(u,v). If v is not a vertex then 
there exists at least one traced ray in B that is incident onto a point v' belonging to edge e such that 
d{u,v') < d(u,v)( 1 + rje'). Here, r] = 2(1 + —) and 0 cm is the maximum value of the critical 

angle at any edge in V. 


Proof: Consider the Type-1 ray r from u to v, which is not traced and a ray r' neighboring r in 
the angular order of rays originated at u is traced by the algorithm. (See Fig. [8|) 

Recall that the angle between any two 
successive rays is 5. We consider the 
weighted Euclidean distance of points on 
the segment (v,v'). If d(u,v') < d(u,v ) 
then we are done. If not, then consider the 

point p on the ray r', at distance d(u,v) 
from u. First, assume that p and v are on Zi? ure 8: Illustrating the construction in proving the Lemma 

EL2] 

a common face. Let w\ be the weight of 

the common face. Thus using Lemma |4.1[ d(p,v ) < d(u,v)e', since the two rays r and r' be¬ 
long to the same bundle. Moreover, the Euclidean distance between p and v is upper bounded 
by d(u,v)f^-. Note that the angle of incidence 6 is less than or equal to the critical angle of e. 





Hence, the vv'p angle is greater than | — 6 cm . Let p' be the projection of p onto vv'. Then 
\\pp'\\ = ||up|| sin pvv 1 = \\pv'\\ sin#'. Thus, \\pv'\\ < and d(v,v') < d(u,v)e '(1 + —). 

Next consider the case in which p and v lie on different faces. Let the line segment vp be crossing 
faces /i, fi ,..., fk. (Here, f\ is the face containing v.) We assume w.l.o.g. that w\ < W 2 ■ ■ ■ < u>k, 
so that rays r and r’ diverge as they progress across these faces. By Corollary [I] the weighted 
Euclidean distance between p and v is upper bounded by d(u,v)e’] and the Euclidean distance 
between p and v is bounded by ■ Furthermore, the line segment vv ' is incident to e at an 

angle larger than the critical angle of e. By similar arguments as above, the bound follows. □ 


Corollary 2 Let u be a vertex and let v be a point on an edge e in V such that the weighted 
Euclidean distance from u to v via a Type-1 ray in bundle B is d(u,v). If v is a vertex then there 
exists at least one traced ray in B that is incident onto a point v 1 belonging to an edge e incident 
to v such that d(u,v') < d(u,v)( 1 + rje'). Here, rj = 2(1 + —) and 0 cm is the maximum value of 

the critical angle at any edge in V. 

The above Lemma and Corollary are used to bound the approximation of the distance from 
u to a point with the condition that the point lies between any two successive rays in the same 
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bundle B. A Type-I weighted shortest path from u to an arbitrary vertex v, however, could use 
a critical source and hence use multiple bundles. The following theorem consider these types of 
paths in ensuring an (1 + e)-approximation as well as the time complexity. 


Lemma 4.3 If the angle between a pair of successive rays is as specified in Lemma \4-l\ then a 
Type-I weighted shortest path from any vertex u to another vertex v in V can be approximated to 
within a factor of (1 + 2 ri 2 ye') using rays in Tr(u), where r/ = 2(1 + 


l 

cos 6 C 


Proof: We first consider when an optimum weighted shortest path from u to v traverses a point p± 
such that a weighted shortest path to p[ using the rays in Tr(u) has a critical source p\ in-between. 
Consider two such successive rays r,; e Tr{u ) and ry 6 Tr(u) that originate from u and that lie in 
the same bundle B (refer Fig. [9]) such that r* is incident to a point p" located on edge e, at an 
angle less than the critical angle of e and rj e Tr{u) is incident to a point p\ located on edge e, at 
an angle equal to the critical angle of e. Due to inequality Q, the condition imposed on e' in the 
assumption of this theorem, and the construction of the rays, such rays are guaranteed to exist in 
the same bundle. 

Let dp(x,y) represent the distance between points x and y along path P. Suppose an optimal 
weighted shortest path P opt intersects edge e between p'[ and p\ at p' x . Let P op t refract at point 
p\ with 6 0 pt as the angle of refraction. Consider a path P that approximates P op t'■ it uses ray rj 
from v\ to pi and then uses a weighted shortest path from p\ to some point V 2 - Then the mod¬ 
ified path has length specified by dp(u,v) < dp(u,pi) + dp(p±,v) where dp(pi,v) is a weighted 


shortest path from p\ to v in the discretized space. By Lemma 4.2, dp{u,p\) < d op t(u,p\)(l + ye') 
with dp(p\.p\) < ye'. Further, d opt (pi,v) < dp(p\,p \) + d op t (p\ ,v). Applying Lemma 4.2 again, 


dp(pi,v) < dopt(pi,v)(l + ye') . Thus dp(u,v) < d op t(v i,^2)(1 + 2 ye'). Since the cardinality of 
any edge sequence of a path is 0(n 2 ), repeating the above analysis for all edges that the ray might 
encounter, results in an error factor of 0{n 2 ye'). The case when no critical source is in-between is 
handled by Lemma 4.2 This completes the proof of Lemma 4.3 □ 



parallel to r op t 


Figure 9: Illustrating the construction in proving Lemma 


4.4 


We next consider Type-II paths. 


Lemma 4.4 If the angle between any pair of successive rays is as specified in Lemma 4-1 an d 
rays generated from a critical segment are e' Euclidean distance apart, then a Type-II weighted 
shortest path from any vertex u to another vertex v in V can be approximated to within a factor of 
(1 + 2 n 2 ye') using rays in Tr(u), where y = 2(1 + —). 


Proof: An optimal Type-II path P can be partitioned into two: a path P\ from v to k and a path 
P 2 from k to w. Let the critical segment k have as a critical point of entry, y on edge e. We first 
show that a good approximation to the path P 2 from n to w can be found. Let P be a weighted 
shortest path that originates at a point p on n and strikes w. This path get critically reflected 
by edge e; hence, exits k. Since rays are generated e 1 apart, there exist a sibling pair of rays that 
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originate from two points p' and p" adjacent to p and strike edges incident to w. These rays are 
parallel to P and a proof similar to Lemma 4.2 shows that P is approximated to well within a 
factor (1 + pe'). Furthermore, the points p' and p" can be discovered by a binary search for the 
point p, with an additive error of e'. Finally, using Lemma [4.3[ it is clear that the weighted shortest 
distance from v to y is approximated to within a factor of (1 + 2 n 2 pe') where d\ is the optimal 
weighted Euclidean distance from v to n. □ 


We finally finish the analysis with the observation that an approximate weighted shortest path 
from s to i can be split into at most n approximate weighted shortest paths. Thus we have the 
following Theorem to summarize. 

Theorem 4.1 Let V be a weighted triangulated polygonal domain with vertex set V. Let e' = 
min{^5^, 4n J m,n i } where p = (1 + —j) and 9 cm is the maximum critical angle of any edge 

inV. A weighted shortest path from sGf to t E V inV can be approximated to within a factor of 
(1 + e) using rays in L) u£ vTr(u) . 


5 Details of the algorithm 

The algorithm is event-driven, where the events are considered by their weighted Euclidean distance 
from s. The algorithm starts with initiating a set TZ(s) of rays that are uniformly distributed around 
s. With each traced ray r, we save the point of its origin, points of refraction, critical points of 
entry and critical points of exit along the path traced by r i.e., as these event points occur in order. 
As r is traced, we append the event points to the list associated with r. The various types of events 
that need to be both determined and handled are described in the following Subsections. 


5.1 Initiating rays from a vertex 


This procedure is invoked to initiate rays from a given vertex, say v, when the discrete wavefront 
strikes v. Since s is also a vertex of V, this procedure is also used in initiating rays from s as well. 

Let f be the face along which a ray has been determined to strike a vertex v on the face. 
Let J-(v) be the collection of all the faces incident to v except for f. The set H(v) of rays are 
initiated from v, and all these rays lie on the set P(v) of faces. Due to Proposition [4] (non-crossing 
property of weighted shortest paths), we do not initiate rays from v over the face f. Further, the 
angle between any two successive rays in TZ(y) at v are bounded by 6 (whose value is bounded as 
described in Section [4]) . 

Consider any face / E -F(u). Let 
vv', vv”, v'v" be the edges bounding /. (See Fig. 

10 ) The set Bf C lZ(v) of rays that lie on face 
/ is a ray bundle: every ray in Bf strikes v'v" 
before striking any other edge in V. For every 
such ray bundle Bf, we find a sibling pair corre¬ 
sponding to Bf. Let r ' (resp. r") be the ray in 
B f that strikes v'v" at g r /(resp. q r „) such that Fi S ure 10: Illustrating successive rays (shown in blue 
,, . _ , .. , ,, color) striking v from /' and a ray bundle being initiated 

there does not exist a ray m Bf that strikes v v c , _ , , „ . ,, .... . c D ,, n 

J J on a iace / £ iF(y); r ,r is the sibling pair oi Bf. When 

between q r / (resp. q r ») and v (resp. v ). We r ' along / and/or r" along f" strike v'v" and v'v"' edges 
do binary search (with respect to edges vv' and respectively, a discrete wavefront is initiated from v'. 
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vv") among the rays in TZ(v) to find both the 

rays r' and r" . We trace the sibling pair r' and r" as long as both the rays refract along the same 
sequence of edges £ until the bundle is forced to split. Let e be the last edge in the edge sequence 
£ at this stage. Let q' and q" be the points on edge e to which rays r' and r" are incident when 
they are traced. Let d v be the weighted Euclidean distance from s to d i.e., when the discrete 
wavefront struck v. The following sets of event points are pushed to the event heap: the event 
point corresponding to tracing the ray r' (resp. r") to q' (resp. q") that occurs at distance equal 
to the weighted Euclidean distance d v added with the weighted distance along r' (resp. r") from 
s. Further, the corresponding sibling pair is saved with each event point. The splitting of sibling 
pairs is further detailed in Subsection |5.4| 

At the initialization step, let /, f" be two adjacent faces in T(v). (See Fig. 10 ) Let vv'v" (resp. 


vv'v'") be the triangle defining / (resp. f"). Let r' be the ray in lZ(y) that lies on / and let r'" 
be the ray in lZ(v) that lies on face f" such that no ray in TZ(v) lies between r' and r'" . Also, let 


q r / (resp. q r m) be the point on edge v'v" (resp. 


v'v'") 


to which ray r' (resp. r'") is incident when 


traced. Then the event point for initiating a discrete wavefront from v' is pushed to the event heap 
with the key value min(d„ + Wf\\vq r '\\ + w v i v "\\q r 'v' \\, d v + Wf \\vq r '" || + w v i v i"\\q r "iv'\\). 

To improve the time complexity, we use the non-crossing property of (weighted) shortest paths 
(Proposition [dj): whenever a ray bundle B that originates at a vertex v' strikes another vertex u, 
we save that information with v so that whenever another ray bundle B' that has originated from 
the same vertex v' splits at vertex v at an event corresponding to a larger distance, we do not split 
B' in order to propagate from v. To achieve this, whenever a vertex v is struck by a ray bundle 
B, we save the origin of B with v in a set. Further, whenever a different ray bundle B' splits due 
to vertex v, we check this set and (i) update the bundle that first strikes v, say B" (ii) eliminate 
progressing any bundle (including B') whose origin is same as B" if it crosses the bundle B". 


5.2 Handling the ray striking an edge critically 

Let r \ , V '2 be a sibling pair and let e(v',v") be a common edge to faces / and f. Consider the 
following event point: A ray r\ is traced along face / and is critically incident to e at point y £ e. 
When this event occurs, we initiate two kinds of rays: set of rays that originate from the 

critical segment yv" (denoted by «:); set 7 Z(y) of Steiner rays that originate from the critical source 
y. In the following Subsections, we describe algorithms to both initiate these sets of rays and to 
set up ray bundles. 


5.2.1 Initiating rays from a critical segment 



Figure 11: Illustrating a sibling pair r,r' originated at Figure 12: Illustrating a sibling pair r,r' originated at 
a critical segment that does not require splitting. a critical segment that does require splitting. 

We describe the procedure to initiate rays from the critical segment n first. Given that a 
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geodesic shortest path can be reflected back onto face / from any point on n, we initiate a discrete 
wavefront from k. Let it be a vector normal to edge e, passing through point y' to some point 
in face f. The initiated rays get critically reflect back into face / from e while making an angle 
—0 c (f, f) with it. Note that the rays in the discrete wavefront originating from k are parallel 
to each other; further, to achieve an (1 + e)-approximation, the distance between any two such 
successive rays along e is upper bounded by e' (here, d is defined as in Lemma 4.3). 

Let 2/1 be the point located on e at Euclidean distance ||yV , || — d from y' such that y\ is located 
between y' and v". Let r and r' be two critically reflected (parallel) rays (making an angle — 0 c (f , f) 
with if), that originate from points y' and y\ respectively. 

If both the rays are incident to the same edge, say e" of /, we set r, r' as a sibling pair and all 
the rays between r and r' that potentially cross the same edge sequence in the future together with 
r,r' is a ray bundle. (See Fig. [Tl}) Let x (resp. x") be the point at which the ray r (resp. r') is 
incident onto edge e" . This sibling pair is traced further across the polygonal domain. 

If the rays in a sibling pair r, r' are incident to distinct edges of /, we need to find two rays 
r\. r' 2 in 'R(k), to respectively pair up with r and r', for forming two new sibling pairs which define 
two corresponding ray bundles. (See Fig. 12 ) Further, we also need to store an event point in the 
event heap that corresponds to the current shortest distance to v if it is via the rays r\ and r' 2 . 
At that stage, a discrete wavefront from v is initiated. The algorithm to split a sibling pair that 
originates from a critical segment is described in Subsection |5.4.1 


5.2.2 Initiating rays from a critical source 

There are three cases to consider based on the origin of the sibling pair, r\ and r2, of a bundle B 
that represents the discrete wavefront from the critical source y’ on an edge e: 

(i) ri,r2 have originated from some critical segment 

(ii) ri,r2 have originated from some vertex 

(iii) ri,r2 have originated from (possibly distinct) nodes of a tree of rays, say Tr(w) 

We first rule out case (i) as it cannot occur, since from Proposition [lj a critical point of exit 
and a critical point of entry cannot occur in succession along a geodesic path. Since case (iii) is a 
generalization of case (ii), herewith we explain event handling for case (iii). 

Let u G Tr(w) be the origin of r\. Note that u may be w itself or a distinct critical source in 
Tr{w). Let r,r' be successive rays in 7 Z(u) such that r is refracted and r' is critically reflected. 
Let y' be the point at which r' is incident to v'v". A node corresponding to critical source y' is 
inserted as a child of u in Tr(vj); further, a set 7 Z{y') of rays are initiated from y'. Let O' be the 
angle at which r refracts (measured with respect to the vector if normal to edge v'v"). Let v'v"v'" 
be the face f onto which ray r is refracted into. Also, let vi and vt, be two vectors in face f" with 
origin y' such that they respectively make 5 and O' angles with respect to if. (See Fig. |4j) The 

rays in TZ(y') are uniformly distributed in the cone C(y ', y'v^' , ut). A ray bundle and a sibling pair 
corresponding to that ray bundle are determined: Let r% be the ray in TZ(y') that subtends the 
minimum required angle, d with edge e. And let r 2 be the ray that originates from y' and is parallel 
to ray r in f" . Then r2,r 3 are extremal pairs of rays in Tr(vj). The bundle B is now modified to 
comprise all the rays that lie between r 2 and r%, with r -2 and r 3 together forming a sibling pair of 
B. This sibling pair is traced over the face f ". 
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When traced, if both and r% are incident to same edge of face f". then they together will 
continue to be a sibling pair. Otherwise, to form two sibling pairs, new rays to pair up with r2 and 


rs are found from the sets of rays initiated in Tr(w ); the procedure is detailed in Subsection 5.4.2 


These sibling pairs’ are traced and the corresponding event points are pushed to the event heap. 


5.3 Extending rays across a region 

Given the bundles of rays that strike an edge e = (a, 6), the rays need to be extended across the 
region they enter. Let T = (a, b, c ) be the triangular region. The rays that refract after striking 
e, as well as the critical rays that originate from e, need to be extended to determine their strike 
points on the other two edges (a, c) and (6, c). 

Extending a bundle may involve splitting a sibling pair. While the split operation details will 
be defined subsequently, we discuss how to manage the extension of the set of bundles. Note that 
only one of the bundles that strike (a, b ), will lead to a weighted shortest path from the origin s to 
the vertex c from amongst the bundles that strike edge (a, b ). Furthermore, from a pair of bundles 
that cross each other, only one will be retained, due to the non-crossing property of shortest paths; 
and, from a pair of bundles that split at vertex c only one bundle need be split, i.e., the one that 
determines the shortest distance to c. 

Determining a weighted shortest path to c via rays contained in a bundle B will be detailed 
later. Given a weighted shortest path to c via each bundle in the current set of bundles B that 
strike e, one can determine a weighted shortest path to c with respect to B. The bundle B' E B 
that determines a weighted shortest path is maintained. Note that B changes as bundles strike 
edge e. The bundle B' partitions the bundles in B into two sets of bundles, one that strike the edge 
(a, c) and the other set that strikes ( b , c). These bundles are traced after processing their strike on 
their corresponding edge. To ensure that the relevant bundles in B are traced across an edge, say 
e = (a, c), we determine a shortest path to the edge e via rays in bundle B E B. When this event 
occurs, the bundle is traced across edge e if it does not violate non-crossing property. A weighted 
shortest path to the edge is determined by a binary search over the space of rays in the bundle. 
This binary search procedure is similar to the process of determining a weighted shortest path to 
a specific vertex by a binary search. 

Note that the above methodology is true for all bundles, independent of the source being a 
vertex or a critical segment. We next determine how to efficiently determine the splitting of sibling 
pairs that define a bundle. 


5.4 Splitting a sibling pair 

When a sibling pair needs to be split i.e., a ray bundle needs to be partitioned into two, then, 
based on the origin of the two rays in the sibling pair being considered, the following procedures 
are invoked. 


5.4.1 Pair that originates from a critical segment 


Consider a sibling pair r, r' that originates from a critical segment n. (See Fig. 12 ) Let £ be the 
edge sequence of r and r' that the two rays have in common till they are incident to an edge e of 
face / and let B be the corresponding ray bundle. Let r be incident to x E e! and let r' be incident 
to x" E e". Given that r and r' are parallel between any two successive edges in £, following are 
the possibilities: r and r' refract from x and x" respectively at the same angle; r and r' critically 
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Figure 13: Illustrating a ray bundle split when that ray bundle is originated from a critical segment k; two new 
sibling pairs are also shown 

reflect from x and x". Since the latter kind of rays are not geodesic paths (Proposition [T]), we focus 
on the former. If the rays refracted from x and x" are both incident to the same edge of /, then 
there is no need to split the sibling pair r,r'. Otherwise, as explained below, we find new rays r[ 
and r 2 that originate from k to respectively pair with r and r 1 . 

Let v \, V 2 , '03 be the vertices of / where the endpoints of e are v \, V 2 - (See Fig. 


13 ) Let l be a 


line segment on face / such that it is parallel to line segment / n r (r refracted at p) and passes 
through vertex V 3 of face /. Let x" G V\V 2 be its other endpoint. We interpolate over rays in 7Z(k) 
to find a ray r[ G 1Z(k) (resp. r' 2 G 7£(ac)) such that the section of r\ (resp. r 2 ) on / lies between l 
and the section of r (resp. r / ) on /. Then the ray bundle B splits into: a bundle with sibling pair 
r,r[ and another with r',r 2 as its sibling pair. Since between any successive edges in £, rays r and 
r' are parallel, an interpolation in possible. Based on the ratio of ||aca: // 1| to ||xa/||, we interpolate to 
find a point q on k so that a critically reflected ray r"' from q reaches V3. The ray r[ (resp. r 2 ) is 
the one whose origin is closest to q among all the rays between r"' and r (resp. r'). These sibling 
pairs’ are pushed to the event heap with the key values being the respective weighted Euclidean 
distances’ from s to the points at which these pairs strike the edges V 1 V 3 and V 2 V 3 . 

Let the ray r[ be incident to V1V3 at z\ and let r 2 be incident to V2V3 at z 2 . An event point to initi¬ 
ate a discrete wavefront from V 3 is pushed to the heap with the key value min(d 2 / + w vlV3 Hz^H, d,/+ 
Wv 2 v 3 II z 2Fill) where d z ^ is the weighted Euclidean distance from s to z\ and d z > is the weighted Eu¬ 
clidean distance from s to z'n. 


5.4.2 Pair that originates from a tree of rays 

Let ej, ej, e*, be the edges bounding face /. Also, let r', r" be a sibling pair of Tr(w ) such that both 
of them strike edge e* of /. With further expansion of the wavefront, let r' strike edge ej at q' and 
let r" strike edge e& at q". This requires us to split the ray bundle, B, corresponding to sibling 
pair r',r". Let P' (resp. P") be the critical ancestor path of r' (resp. r"). Further, let REG be 
the open region bounded by rays along critical ancestor paths P',P", the rays r',r" and the line 
segment q'q". Let R\ (resp. R 2 ) be the set of rays such that a ray r G R\ (resp. r G R 2 ) if and 
only if r originates from a critical point of entry or a critical source located on the critical ancestor 
path P' (resp. P") and the ray r lies in REG. With binary search over the rays in R±, we find a 
ray r\ G R\ and with binary search over R 2 we find a ray r 2 such that r\ intersects e^, V 2 intersects 
efc, and r\ and r 2 are either successive rays originating from the same origin or adjacent origins on 
a critical ancestor path. The binary search is performed over the nodes on the critical path, and 
at each node v the two extreme rays in the set of rays TZ(v) that originate at that node are used 
to decide whether the two rays r\ and r 2 lie within TZ(v) or not. This leads to splitting ray bundle 
B into two B\ and B 2 with r\,r' and r",r 2 sibling pairs, respectively. Events corresponding to 
the ray bundles created are pushed to the event queue reflecting the split. Further, we also push 
the event corresponding to initiating a discrete wavefront from v. Note that splitting a sibling pair 
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that originates from a vertex is just a special case of the procedure listed above. 


6 Improving the time complexity: interpolating versus tracing 
rays 

To improve the time complexity of our algorithm while obtaining a weighted shortest path with 
(1 + e) multiplicative error, instead of tracing a ray r across an edge sequence to find its point (and 
angle) of incidence onto an edge e, we show that it suffices to interpolate the position and the angle 
of refraction of r from the last edge of £. For this, we need an assumption (albeit mild) that the 
angle between any two successive rays in any ray bundle is upper bounded by it — . 


From Lemma 4.1 and Lemma 14.41 we note 


that in order to find an (1 + e)-approximate 
weighted shortest path, the cardinality of the 
set S of rays that originate from a vertex, or a 



critical source, is chosen to be 0( , ). Let 

’ Me')™ +1 

£ be an edge sequence of a sibling pair r / , r". 
and let e be the last edge of £. In doing binary 
search over S to find a ray r G S', we need to 
trace 0(lg(I() + n 2 lg(^-)) rays in S to e. Tak¬ 
ing into account the cardinality of £ (which is 
0(n 2 ) from the Proposition^, the time it takes 
for binary search to find a ray from the source 
is 0(n 2 lg(()) + n 4 lg^)). As detailed below, 
we reduce the time involved in this by interpo¬ 
lation. 

Let r' and r" be two sibling rays in a bundle B in Tr(v), with 4j(= 2r) » +4 ) being the angle 


Figure 14: Illustrating the error in interpolat¬ 
ing versus tracing a ray across an edge sequence 
£ 


between them. (See Fig. 14 ) Consider a ray r £ B with origin v and that lies between r' and 


r". Also, let e p be an edge that is intersected by all three rays r, r', and r", respectively at points 


q, q' and q ". We denote 


M 

I Wq" 


with 7 . (This ratio is known if the target point q on e is known 


even though r itself may not be known.) The interpolation of ray r is denoted with 77 and it is 
characterized as follows: 

(i) Let f p be the face bounded by e p _i and e p . The angle between the vectors induced by f p Cirj 
and f p FI r' is upper bounded by 4(7. 


(ii) The Euclidean distance of interpolated point of incidence qj of ray 77 from q' is q' + 7 ||gV / | 


Lemma 6.1 Let £ be the edge sequence associated with a ray bundle B in Tr(v), having siblings 
r' and r"; and, let e be the last edge in £. Let r G B be a ray that lies in between r' and r" and 
let q be the point of incidence of r on edge e when traced across £. Also, let qj be the point due to 
interpolation. Given that the angle 71 between r' and r" is upper bounded by 2 f n+i , the weighted 
distance d between qj and v is approximated by rj within a multiplicative error of (1 + A), and the 
angle of incidence of rj is (1 + e)-approximated with respect to r. 
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Proof: Let e\, e 2 ,..., e v be the edge sequence 


in £. As in Lemma 4.1 angles of incidence and 
angles of refractions of rays r' and r" are as 
mentioned here. Let 9[, 8 ' 2 ,..., 9' p be the an¬ 
gles at which the ray r' is incident on edges 
e\, e 2 , • • •, e p respectively. (See Fig. 15 ) Let 



Figure 15: Illustrating the angle of incidence 
and refraction of rays r , , r", and r' at an edge 

Si £ £■ 


62 , 63 ,, 6 P+ 1 be the angles at which the ray r' 
refracts at edges e±,e 2 , ■ ■ ■ ,e p respectively. Sim¬ 
ilarly, let 6 \ + d\ , 6’ 2 + d’ 2 ,... , 6 ’ p + dp be the an¬ 
gles at which the ray r" is incident on edges 

e\, e 2 , ■ ■ ■, e p respectively. And, let 62 + 62 , 63 + 63 ,..., 6 p+ \ + 8 p+ \ be the angles at which the ray r" 
refracts at edges e\,e 2 , ■ ■ ■ ,e p respectively. For every i, if a critical angle exists for edge e*, assume 


that both 6 \ and 6 \ + d\ are less than that critical angle. Using Lemma 4.1, we derive a bound on 


5 P as follows: 


< 


< 


TT cos(6>( + (%) 
w p 1 cos(6' i+ i + 5 i+ i) 

wi §r yr cos(fl<)-<%(sin(fl<)) 
w p ' cos 6 ’ i+l 

n ^ 

w p - LJ - 

W\ 


(from Q) 

__i_ 

<5j+i sin 6 i+ i 


cos 6 [ 
COS 6 % 3 ~ 1 


i «si<n 




cos 6 i+ 1 


cos 8 j' 


< 


< 


0 EA)> 


Wp 

^<ns&><n<^ 

^<n +e + o{ ^ «* % - 

* 3 3 


Let 6 [ + 7 i, d ' 2 + 72 , • • ■, Q'p + l'p be the angles at which the ray r is incident on edges e\, e- 2 , ■ ■ ■, e p 
respectively. Also, let 62 + 72 , 83 + 73 ,..., 8 Pl + 7 p+ i be the angles at which the ray r refracts 
from edges ei, e 2 ,..., e p respectively. Then analogous to the above, we can lower bound the angle 
at which ray r strikes edge Bj as follows: 7p > ^ 7 i(n i ^J I )(l - c ^ e , + o(J 2 j cd%))- For 
every 2 < i < p, let fi be the face bounded by ej_i and e t . Also, for every 2 < i < p, let dj be 
the angle between the vectors induced by /,■ n 77 and f ; n r'. The angle dl is upper bounded by 


S-wn, eSStKi + t ., a - 


6i 


jj cos Oj 


+ j cos J e . ). Hence, the ratio of the error in the angle S p versus 
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the angle y p , denoted with R p equals to 


< 


< 


< 


< 


^p~1p 


7 p 


(E 


Si'y 


i cos Oi cos 0; 


- ^V) + °(Ej ^ie~)) 


i + E 


7j 
j cos 9', 




o 2(7 , ^ ^ » o 

n max-— + o( > -— since p < n 

cos Oi ' 


cos #, 


e , 25, 

—r when max 

n 2 


e 

< — 


i cos n 4 


Therefore, to keep the error ratio less than e, we choose d 4 < 2n "+ 2 s i nce min* cos 0,; > (In the 
last inequality, we used the assumption that 0■ < f ~ E)- In fact, we choose 8 \ less than or equal 
to 0 )l n +4 , so that to upper bound R p with e as well as to upper bound \\qi — q ||. 

Also, let q[, q' 2 , ■ ■ ■, Q p = Qi respectively be the points of incidence of interpolated ray 77 on edges 
ei, e 2 ,... , e p in T. Considering that the ray r is traced across the edge sequence of cardinality 0(n 2 ), 
the \\qiq\\ is bounded recursively. Let dist r be the distance along ray r from v to q and let dist ri be 
equal to Ef =2 Iki-i^ll- Since the error in angle accumulates, r J ist is upper bounded by Ay. □ 


In the interpolation method, the binary search is performed on the ordered set of rays, TZ(v), 
at vertex v until the angle between the rays guiding the binary search is less than or equal to „ e +4 . 
This results in 0(n(lg j)) search steps, where each step would require tracing a ray from the source. 
Moreover, this method can also be used to approximate the shortest distance to an edge via rays 
in a ray bundle. We summarize the discussion in the following Lemma. 

Theorem 6.1 The interpolation method determines the shortest distance to a given point p from 
a source v to within a multiplicative approximation factor of (1 + e) in 0 (n lg(")) time using rays 
from a ray bundle. Furthermore, an (1 + e) -approximation to the shortest distance from v to an 
edge e can be determined in 0(n(lg ™)(lg ^)) time using rays from a ray bundle. 


7 Analysis 

The analysis is based on showing the following facts: (i) bundles of rays are correctly maintained, 
and (ii) the shortest distances to vertices, and edges similarly, are correctly computed using the 
bundles. 

Bundles are initiated from vertices and propagated across the faces of the domain as they strike 
the edges of the domain. The propagation is evidently correct, the modifications to the bundles 
being (a) bundle splits and (b) elimination of bundles due to crossing of bundles. As described in 
Subsection |5.3t elimination of bundles is determined at vertices when more than one of the bundles 
split at that vertex. Note that bundles are propagated until they split and further propagation 
of the split parts of the bundles occurs when the shortest distance event, of weighted distance d, 
corresponding to the bundle, say B, striking an edge is determined by the heap. When two bundles 
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cross each other completely, one of them is eliminated from further consideration. Thus for the rest 
of the proof, we will assume w.l.o.g., that bundles are correctly maintained. 

To determine approximate weighted shortest paths to vertices correctly, we consider two cat¬ 
egories of paths, Type-I and Type-II. Recall that Type-I paths do not have critically reflected 
segments in between, while Type-II paths does have. Analogous arguments prove the correctness 
involved in computing approximate weighted shortest paths to edges. 

Lemma 7.1 The algorithm correctly determines sibling pair in Tr(u) for every source u, and 
computes an (1 + e)-approximation to a Type-I weighted shortest path from u to v, for any vertex 
v in V. 


Proof: Suppose that there exists a Type-I weighted shortest path between u and v, traversing an 
edge sequence £. Lemma 4.3 shows that an (1 + e) approximation can be found using rays within a 
bundle. Bundles are maintained using sibling pair. We show that the algorithm maintains sibling 
pair for edges in £ that determine bundles, the sibling pairs being traced rays such that these rays 
can be used to find an (1 + e)-approximate weighted shortest path from u to v. 

Initially when rays are generated from a vertex, v, siblings are computed correctly and each 
ray is refracted correctly. Consider the procedure that computes a sibling pair r\ , r 2 closest to a 
vertex v. The rays r\ and ?’2 are thus two successive rays and the distance to v via rays r± and 
r 2 is updated. Suppose both the rays r\ and r 2 have the same origin. Using binary search with 
interpolation, we find the siblings correctly. Otherwise, the origin of n is different from the origin 
of r 2 - In this case, we find a sibling pair via the binary search with interpolation on the rays in 
the critical ancestor paths of two siblings, say r' and r" such that the rays r\ and 7*2 are part of 
the set of rays lying in between r' and r". The correctness of the interpolation method is shown in 
Lemma 16.11 □ 


Lemma 7.2 Let P be a Type-II weighted shortest path from a vertex v to another vertex w on V 
with a critical segment k in-between. Then the algorithm determines a pair of traced rays in 7Z(v) 
that can approximate a weighted shortest path, P. 


Proof: An optimal path P can be partitioned into sub-paths, each such sub-path going from one 
vertex v to another vertex w. Suppose a sub-path uses a critical segment tz and is partitioned as 
follows: a path P\ from v to k and a path P 2 from k to w. Let e be the edge on a face / such 
that the critical segment k lies on e and reflects rays back onto face /. Also, let y be the critical 
point of entry into k. The correct determination of the critical point of incidence from a bundle of 
rays follows from Lemma 4.3 If w is the endpoint of k then we are done since the distance to the 
endpoint from the critical source is included in consideration. Otherwise rays are generated from 
k that are parallel and separated by a small weighted Euclidean distance less than e'. Let r 1 , 7-2 be 
sibling rays that originate from k such that w lies in between r\ and r 2 - Lemma |4.4| shows that the 
shortest distance from k to w can be found by tracing the rays r\ and 7-2 and interpolating between 
them to find the point on k closest to w. □ 


Theorem 7.1 The algorithm computes an (1 + e)-approximate weighted shortest path from s to t 
in 0(n 5 ( lg”)(lg^)) time. 
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Proof: The correctness follows from Lemmas 17.11 and 17.21 

We first bound the number of ray bundles. The number of bundles that are initiated from 
vertex sources are 0(n). By Proposition 3, the number of bundles that are initiated from critical 
segments is bounded by 0(n 2 ). The processing of such bundles will be handled separately. 

The number of bundles increase when split. The increase is charged to the vertex v that causes 
the split. Only one bundle, the bundle that determines an approximate weighted shortest distance 
to v is split and thus the total number of bundles split are 0(n). Bundles are propagated further 
into adjacent faces after striking an edge. 

We thus need to analyze the time involved in finding an approximate weighted shortest paths 
and splitting ray bundles. Let S be the set of rays from a vertex source v\. To search for a ray 
in S that is closest to the vertex or provides the shortest distance to the last edge e of an edge 
sequence £, the algorithm utilizes the combination of binary search and interpolation (as described 
in Section [6]). Using Theorem 6.1, the complexity of determining the shortest distance when rays 
from a specific vertex are considered is equal to 0{n 3 lg(f) lg(^)). Note that 0(n 2 ) edges will be 
encountered as rays are traced. 

Next, let us consider the case in which the sibling pair is not from a vertex but arises from a 
critical source contained in a tree of rays, say Tr{u). Let S be the set consisting of all the critical 
points of entry in any critical ancestor path P in Tr(u). A binary search on the critical path 
determines vertex u such that the ray to be determined is in 7 Z(u). Since S is 0(n 2 ) in size and 
since it is ensured that no two ray bundles cross, binary search can be used to find a vertex v in 
Tr(u ) such that an approximate weighted shortest path lies in TZ(v) with 0(lg(n)) ray tracings. A 
search on the rays from that specific vertex in 1Z(u) follows as described in the earlier paragraph. 

A similar search is required to determine the ray that identifies the approximate weighted 
shortest distance to an edge, which is required when a bundle splits. The procedure and time 
complexity of this is similar to the determination of a weighted shortest path to a vertex. 

Since the total number of vertices and critical sources is 0(n 2 ), the total work in splitting and 
initiating ray bundles from these sources is 0(n 5 lg(^) lg(4)=)), as claimed in Theorem 


6.1 


Determining the ray originating from a critical segment and that strikes a vertex can be deter¬ 
mined via interpolation (or a binary search) on the space of the parallel rays that are part of the 
bundle that originates from a critical segment. The work involved in tracing a ray that originates 
at a critical segment and traverses an edge sequence takes 0(n 2 ) time (Proposition^. The specific 
pair of successive rays of interest can be found by interpolation, taking 0(1) time. As there are 
0{n) vertices and 0(n 2 ) critical segments (Proposition 3), the time complexity is 0(n 5 lg(^)), 
including the binary search involved. □ 


Single-source approximate weighted shortest path queries 

Here we devise an algorithm to preprocess V to construct a shortest path map so that for any given 
query point q in V an approximate weighted shortest path from s can be computed efficiently. 
As part of preprocessing, using the non-crossing property of shortest paths (Proposition [4]), we 
compute for every edge e in the subdivision a minimum cardinality set 13(e) of bundles such that 
for every point p E e an approximate weighted shortest path can be computed from a bundle 
B E 13(e). Further, for every edge e, the set 13(e) of bundles are saved with e. For a query point 
q E V, exploiting Lemma 4.2 and Corollary [2j we locate q in the triangulation and determine an 
approximate weighted shortest path from q to s via one of the bundles associated to edges of the 
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triangle containing q. 

As mentioned, let N be the maximum coordinate value used in describing V. Since the number 
of bundles is 0 (n 2 ) and since it takes 0(n 2 (lg f )(lg ^)(lgiV)) time to determine the bundles that 
need to be associated to an edge and there 0(n) edges, it requires 0(n 5 (lg ™)(lg ^=)(lgiV)) time for 
preprocessing. The query phase takes 0(n 4 (lg y)(lg ^)(lgiV)) time: there could be 0(n 2 ) bundles 
associated to edges of the triangle containing q and finding an approximate weighted shortest path 
via any one such bundle takes 0(n 2 (lg y)(lg ^)(lgiV)) time. 


8 Conclusions 

In this paper we have presented a polynomial-time algorithm for finding an approximate weighted 
shortest path between two given points s and t. The main ideas of this algorithm rely on progressing 
the discretized wavefront from s to t. The time complexity of our algorithm is 0(n 4 (nlgn + 
polylog (n,e,9mi n ))). Significantly, our algorithm is polynomial with respect to input parameters. 
This result is about a cubic factor (in n) improvement over the Mitchell and Papadimitriou’s 
’91 result [8J in finding a weighted shortest path between two given points, which is the only 
known polynomial time algorithm for this problem to date. In addition, we extend our algorithm 
to answer single-source weighted shortest path queries. Further, with minor modifications, our 
algorithm appears extendable to determine geodesic shortest paths on the surface of a 2-manifold 
whose faces are associated with positive weights. Since the number of events in the problem stand 
at II(n 4 ) (from [8]), it would be interesting to explore further improvements in devising a more 
efficient polynomial time approximation scheme. 
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